Matlab Codes

Chapter 4

Matlab code 4.2: Matlab file “Figure 4-4.m”

%--------------------------------------------------------------------

% This code can be used to generate Figure 4.4

%--------------------------------------------------------------------

 

clear all;

close all;

%% the micro-Doppler signature theoretical curve of the rotating target(narrow band bistatic radar)

c = 3e8;

j = sqrt(-1);

fc = 10e9; % carrier frequency of transmitted signal

v = 0; % translational velocity of target

radt = [0 0 0]; % coordinates of transmitter

radr = [7000 0 0]; % coordinates of receiver

cord = 1000*[3 4 5]; % coordinates of target centre in the world coordinate system

colo = [0 0 0]; % local coordinate system's origin in the world coordinate system

w = [pi 2*pi pi]; % angular velocity of rotation

r = [0.2 0.2 -0.6;-0.5950 0.1828 0.2293;0.3950 -0.3828 0.3707]; % rotation radius vector of the scaterers

ae = [0 pi/4 pi/5]; % Euler angle

ri1 = [cos(ae(1)) -sin(ae(1)) 0;sin(ae(1)) cos(ae(1)) 0;0 0 1];

ri2 = [1 0 0;0 cos(ae(2)) -sin(ae(2));0 sin(ae(2)) cos(ae(2))];

ri3 = [cos(ae(3)) -sin(ae(3)) 0;sin(ae(3)) cos(ae(3)) 0;0 0 1];

ri = ri1*ri2*ri3; % initial rotation matrix

t = 1; % radar illumimated time

prf = 1000; % pulse repetition frequency

pri = 1/prf; % pulse repetition interval

dt = 0:pri:t-pri; % time sampling interval

m = length(dt);

nb0 = (radt-cord)/sqrt(sum((radt-cord).^2))+(radr-cord)/sqrt(sum((radr-cord).^2));

tn = length(r(1,:)); % number of scatterers

r0 = zeros(3); % scatterers in the local coordinate system

for i = 1:tn

    r0(i,:) = r(i,:)-colo;

end

r0r = ri*r0'; % scatterers in the reference coordinate system

w = ri*w'; % angular velocity of rotation in the reference coordinate system

omega = sqrt(sum(w.^2)); % angular frequency

we = w/omega; % unit vector of rotation

wr = [0 -we(3) we(2);we(3) 0 -we(1);-we(2) we(1) 0]; % skew symmetric matrix

fmd = zeros(tn,m); % theoretical micro-Doppler frequency

for i = 1:m

    fmd(:,i) = (omega*fc/c)*(wr*(wr*sin(omega*dt(i))+eye(3)*cos(omega*dt(i)))*r0r)'*nb0';

end

 

%%  the time-requency analysis result of the rotating target micro-Doppler characteristic(narrow band bistatic radar)

rt = zeros(tn,length(dt));% distance between the scatterers and radar

for i = 1:m

    for n = 1:tn

        rt(n,i) = sqrt(((radt-cord)'+v*t+((eye(3)+wr*sin(omega*dt(i))+wr^2*(1-cos(omega*dt(i))))*r0r(:,n)))'...

                      *((radt-cord)'+v*t+((eye(3)+wr*sin(omega*dt(i))+wr^2*(1-cos(omega*dt(i))))*r0r(:,n))))...

                 +sqrt(((radr-cord)'+v*t+((eye(3)+wr*sin(omega*dt(i))+wr^2*(1-cos(omega*dt(i))))*r0r(:,n)))'...

                      *((radr-cord)'+v*t+((eye(3)+wr*sin(omega*dt(i))+wr^2*(1-cos(omega*dt(i))))*r0r(:,n))));

    end

end

s = zeros(length(dt),tn); % echo signal matrix

for i = 1:tn

    s(:,i) = exp(j*2*pi*fc*2*rt(i,:)'/c);

end

st = sum(s,2);

N = 200; % number of Gabor coefficients in time

Q = 50; % degree of oversampling

tfr = tfrgabor(st.',N,Q); % Gabor transform

tfr_r = fftshift(tfr,1);

figure

contour(linspace(min(dt),max(dt),length(tfr_r(1,:))),linspace(-prf/2,prf/2-1,length(tfr_r(:,1))),tfr_r,50)

xlabel('Time (s)')

ylabel('Frequency (Hz)')

axis([0,1,-300,300])